convert p values to qvalues

#read in all phase1 significance values:
all.sigs<-read.csv(file="~/Desktop/anoph.3.march.2020/all.sigs.csv")
table(all.sigs$chrom)
## 
##     2L     2R     3L     3R      X 
## 638059 946059 524268 737401 198852
#convert 0 bayescenv q-vals to 1e-4, to avoid infinite values when plotting with -log10()
all.sigs$prec.q[all.sigs$prec.q == 0]<-.0001
all.sigs$temp.q[all.sigs$temp.q == 0]<-.0001
all.sigs$hum.q[all.sigs$hum.q == 0]<-.0001

#convert lfmm prec p values to genome-wide FDR adjusted q values
for (i in levels(all.sigs$chrom)){
qvals<-qvalue(all.sigs$prec.p[all.sigs$chrom== i])
all.sigs$prec.p[all.sigs$chrom== i]<-qvals[["qvalues"]]
}
#convert lfmm hum p values to genome-wide FDR adjusted q values
for (i in levels(all.sigs$chrom)){
  qvals<-qvalue(all.sigs$hum.p[all.sigs$chrom== i])
  all.sigs$hum.p[all.sigs$chrom== i]<-qvals[["qvalues"]]
}
#convert lfmm temp p values to genome-wide FDR adjusted q values
for (i in levels(all.sigs$chrom)){
  qvals<-qvalue(all.sigs$temp.p[all.sigs$chrom== i])
  all.sigs$temp.p[all.sigs$chrom== i]<-qvals[["qvalues"]]
}
rm(qvals) #space saver

#give order for plotting manhattan plots
all.sigs$order<-c(all.sigs$pos[all.sigs$chrom== "2R"],
                  all.sigs$pos[all.sigs$chrom== "2L"]+61545105,
                  all.sigs$pos[all.sigs$chrom== "3R"]+61545105+49364325,
                  all.sigs$pos[all.sigs$chrom== "3L"]+61545105+49364325+53200684,
                  all.sigs$pos[all.sigs$chrom== "X"]+61545105+49364325+53200684+41963435)

lfmm manhattans

#LFMM precip plot -log10(signficance vals) as manhattan plots
plot1<-ggplot(data=all.sigs, aes(x=order, y=-log10(prec.p), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant precip outliers
lfmm.prec.outliers<- all.sigs[all.sigs$prec.p < .01,]
nrow(lfmm.prec.outliers) #check number of outliers
## [1] 17113
#LFMM hum plot -log10(signficance vals) as manhattan plots
plot2<-ggplot(data=all.sigs, aes(x=order, y=-log10(hum.p), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant precip outliers
lfmm.hum.outliers<- all.sigs[all.sigs$hum.p < .01,]
nrow(lfmm.hum.outliers) #check number of outliers
## [1] 10981
#LFMM temp plot -log10(signficance vals) as manhattan plots
plot3<-ggplot(data=all.sigs, aes(x=order, y=-log10(temp.p), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant precip outliers
lfmm.temp.outliers<- all.sigs[all.sigs$temp.p < .01,]
nrow(lfmm.temp.outliers) #check number of outliers
## [1] 12396
#
grid.arrange(plot1,plot2,plot3, nrow=3)

#bayescenv manhattans

#bayescenv precip plot -log10(signficance vals) as manhattan plots
plot1<-ggplot(data=all.sigs, aes(x=order, y=-log10(prec.q), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant precip outliers
baye.prec.outliers<- all.sigs[all.sigs$prec.q < .01,]
nrow(baye.prec.outliers) #check number of outliers
## [1] 3857
#bayescenv hum plot -log10(signficance vals) as manhattan plots
plot2<-ggplot(data=all.sigs, aes(x=order, y=-log10(hum.q), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant hum outliers
baye.hum.outliers<- all.sigs[all.sigs$hum.q < .01,]
nrow(baye.hum.outliers) #check number of outliers
## [1] 3688
#bayescenv temp plot -log10(signficance vals) as manhattan plots
plot3<-ggplot(data=all.sigs, aes(x=order, y=-log10(temp.q), col = chrom))+
  geom_scattermore(pointsize = 1.2)+
  theme(panel.grid=element_blank(), panel.background=element_blank())+
  geom_hline(yintercept=-log10(.01), linetype="dashed")
#pull all significant temp outliers
baye.temp.outliers<- all.sigs[all.sigs$temp.q < .01,]
nrow(baye.temp.outliers) #check number of outliers
## [1] 10195
grid.arrange(plot1,plot2,plot3, nrow=3)

par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(na.omit(lfmm.hum.corrs$id[abs(lfmm.hum.corrs$diff) < .05]))
sig.pos<-sig.pos[1:27]
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%
par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(na.omit(lfmm.hum.corrs$id[abs(lfmm.hum.corrs$diff) > .8]))
sig.pos<-sig.pos[1:27]
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%

check out the correlations of allele frequency patterns with humidity

dev.off()
## null device 
##           1
#visualize
plot(lfmm.hum.corrs$diff)
abline(h=0, col="red")
plot(lfmm.hum.corrs$phase1, lfmm.hum.corrs$phase2)
hist(lfmm.hum.corrs$diff)
plot(-log10(lfmm.hum.outliers$hum.p), all.cases$diff)
plot(-log10(lfmm.hum.outliers$hum.q), all.cases$diff)
plot(lfmm.hum.outliers$hum.fst, all.cases$diff)

different classes of relationships?

dev.off()
## null device 
##           1
layout(matrix(c(1,1,1,1,2,3,4,5),2, 4, byrow=TRUE))
plot(lfmm.hum.corrs$phase1)
abline(h=c(.4,.1,-.25), col="red")
hist(lfmm.hum.corrs$diff[lfmm.hum.corrs$phase1 > .4])
abline(v=c(-.2,.2), col="red")
hist(lfmm.hum.corrs$diff[lfmm.hum.corrs$phase1 < .4 & lfmm.hum.corrs$phase1 > .1])
abline(v=c(-.2,.2), col="red")
hist(lfmm.hum.corrs$diff[lfmm.hum.corrs$phase1 < .1 & lfmm.hum.corrs$phase1 > -.25])
abline(v=c(-.2,.2), col="red")
hist(lfmm.hum.corrs$diff[lfmm.hum.corrs$phase1 < -.25])
abline(v=c(-.2,.2), col="red")

strongly positive correlation

par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(lfmm.hum.corrs$id[lfmm.hum.corrs$phase1 > .4])
sig.pos<-sample(sig.pos, size=27, replace = FALSE)
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%

weakly positive correlation

par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(lfmm.hum.corrs$id[lfmm.hum.corrs$phase1 < .4 & lfmm.hum.corrs$phase1 > .1])
sig.pos<-sample(sig.pos, size=27, replace = FALSE)
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%

weakly negative correlation

par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(lfmm.hum.corrs$id[lfmm.hum.corrs$phase1 < .1 & lfmm.hum.corrs$phase1 > -.25])
sig.pos<-sample(sig.pos, size=27, replace = FALSE)
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%

strongly negative correlation

par(mfrow=c(3,3))
#plot
sig.pos<-as.vector(lfmm.hum.corrs$id[lfmm.hum.corrs$phase1 < -.25])
sig.pos<-sample(sig.pos, size=27, replace = FALSE)
j<-1 #initialize loop tracker
pb <- txtProgressBar(min = 0, max = length(sig.pos), style = 3) #initialize progress bar
## 
  |                                                                            
  |                                                                      |   0%
for (i in sig.pos){
  frq<-as.numeric(as.vector(phase1.all.freqs[phase1.all.freqs$id == i,3:16]))
  frq2<-as.numeric(as.vector(phase2.all.freqs[phase2.all.freqs$id == i,3:23]))
  plot(sort.pops.phase1$hannual,frq, main = i, ylim =c(0,1))
  abline(lm(frq~sort.pops.phase1$hannual))
  points(sort.pops.phase2$hannual,frq2, col="red")
  abline(lm(frq2~sort.pops.phase2$hannual), col="red")
  setTxtProgressBar(pb, j)
  j<-j+1
}
## 
  |                                                                            
  |===                                                                   |   4%
## 
  |                                                                            
  |=====                                                                 |   7%
## 
  |                                                                            
  |========                                                              |  11%
## 
  |                                                                            
  |==========                                                            |  15%
## 
  |                                                                            
  |=============                                                         |  19%
## 
  |                                                                            
  |================                                                      |  22%
## 
  |                                                                            
  |==================                                                    |  26%
## 
  |                                                                            
  |=====================                                                 |  30%

## 
  |                                                                            
  |=======================                                               |  33%
## 
  |                                                                            
  |==========================                                            |  37%
## 
  |                                                                            
  |=============================                                         |  41%
## 
  |                                                                            
  |===============================                                       |  44%
## 
  |                                                                            
  |==================================                                    |  48%
## 
  |                                                                            
  |====================================                                  |  52%
## 
  |                                                                            
  |=======================================                               |  56%
## 
  |                                                                            
  |=========================================                             |  59%
## 
  |                                                                            
  |============================================                          |  63%

## 
  |                                                                            
  |===============================================                       |  67%
## 
  |                                                                            
  |=================================================                     |  70%
## 
  |                                                                            
  |====================================================                  |  74%
## 
  |                                                                            
  |======================================================                |  78%
## 
  |                                                                            
  |=========================================================             |  81%
## 
  |                                                                            
  |============================================================          |  85%
## 
  |                                                                            
  |==============================================================        |  89%
## 
  |                                                                            
  |=================================================================     |  93%
## 
  |                                                                            
  |===================================================================   |  96%

## 
  |                                                                            
  |======================================================================| 100%